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ABSTRACT 

We present a detailed analysis of our ability to determine stellar masses, ages, reddening and extinction 
values, and star formation rates of high-redshift galaxies by modeling broad-band SEDs with stellar 
population synthesis. In order to do so, we computed synthetic optical-to-NIR SEDs for model galaxies 
taken from hydrodynamical merger simulations placed at redshifts 1.5 < z < 2.9. Viewed under 
different angles and during different evolutionary phases, the simulations represent a wide variety of 
galaxy types (disks, mergers, spheroids). We show that simulated galaxies span a wide range in SEDs 
and color, comparable to these of observed galaxies. In all star- forming phases, dust attenuation has 
a large effect on colors, SEDs, and fluxes. The broad-band SEDs were then fed to a standard SED 
modeling procedure and resulting stellar population parameters were compared to their true values. 
Disk galaxies generally show a decent median correspondence between the true and estimated mass 
and age, but suffer from large uncertainties (AlogM = — O.OGto^ia, A\og age^ — +0 .OSt^o'^) . During 



the merger itself, we find larger offsets: AlogM = -0.13I^;i>nd A logage^ = -0.121^;^^. E{B-V) 
values are generally recovered well, but the estimated total visual absorption Ay is consistently too 
low, increasingly so for larger optical depths {AAy = — 0.54lQ4g in the merger regime). Since the 
largest optical depths occur during the phases of most intense star formation, it is for the highest 
SFRs that we find the largest underestimates (A log S'i^i? — — 0.44j^Q'3^ in the merger regime). The 
masses, ages, E{B — V), Ay, and SFR of merger remnants (spheroids) are very well reproduced. 
We discuss possible biases in SED modeling results caused by mismatch between the true and template 
star formation history, dust distribution, metallicity variations and AGN contribution. Mismatch 
between the real and template star formation history, as is the case during the merging event, drives 
the age, and consequently mass estimate, down with respect to the true age and mass. However, the 
larger optical depth toward young stars during this phase reduces the effect considerably. Finally, we 
tested the photometric redshift code EAZY on the simulated galaxies placed at high redshift. We find 
a small scatter in A2:/(l -|- z) of 0.031. 

Subject headings: galaxies: distances and redshifts - galaxies: high redshift - galaxies: ISM - galaxies: 
stellar content 
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1. INTRODUCTION 

Understanding the growth and aging of galaxies over 
cosmic time requires reliable estimates of their mass, for- 
mation epoch and star formation history. With the cur- 
rent generation of telescopes, stellar velocity dispersion 
measurements can probe the gravitational potential in 
which the baryonic galaxy content resides out to z ^ 1.3 
(van Dokkum & Stanford 2003; Holden et al. 2005). Be- 
yond this redshift, gas velocity dispersions can be mea- 
sured from emission lines, but do not always trace the 
potential due to outflows (Franx et al. 1997; Pettini ct 
al. 1998, 2001; Shapley et al. 2003), and would lead to 
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biased samples missing quiescent galaxies lacking emis- 
sion lines in their spectra (Kriek et al. 2006). For these 
reasons, most studies of high-redshift galaxies have used 
stellar mass estimates derived by modeling of the broad- 
band stellar energy distribution (SED) to characterize 
the mass. 

Since age estimates from Ha equivalent widths (van 
Dokkum et al. 2004; Erb et al. 2006c) or Balmer/4000A 
break strengths (Kriek et al. 2006) are very demand- 
ing in terms of telescope time and only attainable for 
the brightest galaxies, stellar ages as well are commonly 
derived from broad-band photometry. 

Over the past few years, SED modeling has been 
proven extremely valuable in characterizing the galaxy 
population in the early universe (e.g., Papovich et al. 
2001; Shapley et al. 2001, 2005; Forster Schreiber et 
al. 2004). Nevertheless, a number of assumptions are 
required for the limited number of data points (11 pass- 
bands in our case, but often less) to lead to a single solu- 
tion in terms of physical properties such as stellar mass, 
stellar age, dust extinction, star formation rate (SFR), 
and often redshift. 

First, the star formation history (SFH) is generally 
modelled by a simple functional form: a single burst, 
constant star formation, or an exponentially declining 
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model. In reality, liigh-redsliift galaxies show evidence of 
more complex SFHs, often with brief recurrent episodes 
of star formation (e.g. Papovich et al. 2001; Fergu- 
son et al. 2002; Papovich et al. 2005). Second, we 
use the approximation of a single foreground screen of 
dust in accounting for the attenuation, even though in 
reality the dust will be distributed in between the stars. 
Third, we fit solar metallicity models. This is a com- 
mon approach in modeling of high-redshift galaxy SEDs 
in order to reduce the number of degrees of freedom in 
the fitting procedure by one (e.g., Forster Schreiber et 
al. 2004; Shapley et al. 2005; Rigopoulou et al. 2006). 
In addition, the tracks and spectral libraries on which 
the stellar population synthesis models are based, have 
the best empirical calibration for solar metallicity. Al- 
though consistent with the current metallicity estimates 
from near-infrared (NIR) spectroscopy of high-redshift 
galaxies (van Dokkum et al. 2004; Erb et al. 2006a; 
Maiolino et al. 2008), it must be kept in mind that these 
measurements are currently limited to the bright end of 
the galaxy population. Fourth, SED modeling generally 
assumes a piirely stellar origin of the light, while obser- 
vational evidence for a substantial fraction of low lumi- 
nosity AGN at high redshift has been accumulating (van 
Dokkum et al. 2004; Reddy et al. 2005; Papovich et al. 
2006; Kriek et al. 2007; Daddi et al. 2007b). They may 
contribute to the optical SEDs. 

Finally, one adopts a certain attenuation law, initial 
mass function (IMF), and stellar population synthesis 
code. Their appropriateness at low and high redshifts is 
much debated. 

In this paper, we address the impact of the first four 
assumptions (related to SFH, dust attenuation, metal- 
licity, and AGN) using hydro dynamical simulations of 
merging galaxies (see Robertson et al. 2006; Cox et al. 
2006). The SPH simulations follow the star formation on 
a physical basis, resulting in more complex SFHs than 
are allowed in typical SED modeling. They keep track of 
the distribution and metallicity of gas and stellar parti- 
cles, allowing a determination of the line-of-sight depen- 
dent extinction toward each stellar particle separately 
and a knowledge of the stellar metallicity as a function 
of time. Here, we apply the same SED modeling that we 
use for observed galaxies to broad-band photometry ex- 
tracted from the simulation outputs, and study how well 
the mass, age, dust content, and SFR of the simulated 
galaxies can be recovered. 

The reason we use merger simulations for this exer- 
cise is threefold. First, galaxy mergers tiiQ believed to 
play an important role in galaxy evolution (see, e.g., 
Holmberg 1941; Zwicky 1956; Toomre & Toomre 1972; 
Toomre 1977), increasingly so at high redshift (see, e.g., 
Glazebrook et al. 1995; Driver, Windhorst, & Griffiths 
1995; Abraham et al. 1996). Moreover, along their evo- 
lutionary path they are visible as vastly different galaxy 
types, allowing to test the recovery of stellar population 
parameters under a wide range of conditions: gas-rich 
star-forming disks, dust-obscured mergers, and quiescent 
spheroids. Finally, in a separate paper we will compare 
predictions of the color distribution and mass density 
of high-redshift galaxies derived from these simulations 
with the observed galaxy population in deep fields. A 
good understanding of what it is we measure with SED 
modeling is crucial in order to compare identical mass- 



limited samples of observed and simulated galaxies. 

We start with a description of the simulations in §2. 
Next, we explain the methodology of our SED modeling 
in §3. §4 discusses how well we can measure stellar pop- 
ulation properties when a spectroscopic redshift is avail- 
able. §5 repeats the analysis, now leaving the redshift as 
an extra free parameter (i.e., fitting for the photometric 
redshift). Finally, we summarize the results in §6. 

2. THE SIMULATIONS 
2.1. Main characteristics 

The simulations on which we test our SED modeling 
were performed by Robertson et al. (2006). We refer 
the reader to that paper for a detailed description of 
the simulations. Briefly, the simulations were performed 
with the parallel TreeSPH code GADGET-2 (Springel 
2005). The code uses an entropy-conserving formulation 
of smoothed particle hydrodynamics (Springel & Hern- 
quist 2002), and includes gas cooling, a multiphase model 
for the interstellar medium (ISM) to describe star for- 
mation and supernova feedback (Springel & Hernquist 
2003), and a prescription for supermassive black hole 
growth and feedback (Springel et al. 2005b). 

At the start, each simulation consists of 120000 dark 
matter particles, 80000 gas particles, and 80000 stellar 
particles. They represent two stable, coplanar disk galax- 
ies, each embedded in an extended dark matter halo with 
Hernquist (1990) profile. We have realisations where the 
disks start with a gas fraction of 40% and 80%. Stel- 
lar masses at the start of the simiilation varied from 
7.0 X 10^ Mq to 2.3 X 10" Mq per disk galaxy. In total, 
we study 6 simulations (2 different gas fractions, and 3 
different masses) for which at 30 timesteps all informa- 
tion was stored in a snapshot. Each snapshot was ana- 
lyzed as seen from 30 different viewing angles and placed 
at 8 redshifts between z = 1.5 and 3 (§2.2). We repeated 
our analysis on a limited number of snapshots and view- 
ing angles of merger simulations of non-coplanar disks 
and conclude that the results presented in this paper are 
robust against such variations. 

For a given virial velocity, the halo concentration, virial 
mass and virial radius were scaled following Robertson et 
al. (2006) to approximate the structure of disk galaxies 
at redshift z = 3. In practice, this means that the mass- 
and redshift-dependent halo concentration measured by 
Bullock et al. (2001) was adopted: 



M, 



collfi 



-0.13 



(1) 



where McoU,o ~ 8 x 10 Mq is the linear collapse 

mass at z=0, and that the following scaling relations were 
used for the virial mass and virial radius of the progeni- 
tors: 
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where Vyir is the virial velocity and H{z) is the Hubble 
parameter. Disk sizes were initialized according to the 
Mo et al. (1998) formalism for dissipational disk galaxy 
formation assuming the fraction of the total angular mo- 
mentum contained in the disk equals the fraction of the 



3 



total mass contained in the disk ( 



0.041 to match 



the Milky Way-like model by Springel et al. 2005a). The 
disk scale length is then derived from the halo concen- 
tration Cvir (Eq. 1) and the galaxy spin A. The adopted 
value of A = 0.033 is motivated by cosmological N-body 
simulations (Vitvitska et al. 2002). 

We set the ages of the stars existing at the start of 
the simulation such as to represent a constant star for- 
mation history prior to the start of the simulation at 
a star formation rate (SFR) equal to that calculated in 
the first phases of the simulation. The corresponding 
metallicities of stars present at the start of the simula- 
tion were then set according to the closed box model: 
Z{t) = —yln[fgas{t)], where Z{t) is the metallicity of a 
stellar particle formed at time t, the yield y=0.02 and 
fgas{t) is the gas fraction of the system at the considered 
time. Similarly, the gas at the start of the simulation was 
assigned a uniform metallicity Zgas{ts) = -yM.fga.s{ts)] 
where tg represents the start of the simulation, and 
,fgas{ts) = 0.4 or 0.8 respectively for our 2 gas fraction 
runs. The closed box model represents an upper limit 
on the allowed enrichment by heavy elements, which in 
reality may be reduced by outflows or infall of metal- 
poor gas (Edmunds 1990). The fact that we consider 
2 gas fractions guarantees a wide range of progenitor 



types, with ages of a few 100 Myr and Zg 



0.004 



at the start of the simulation for fgas = 0.8 to typical 
stellar ages of a Gyr and nearly solar gas metallicity for 
fgas = 0.4. Starting from the initial conditions described 
above, the GADGET-2 code computed the subsequent 
evolution of the stellar populations using the star for- 
mation and metal enrichment prescriptions outlined by 
Springel & Hernquist (2003). The mass- weighted mean 
stellar metallicity for the fgas = 0.4 runs increases in 
- 700 Myr from 0.4 Zq to 0.9 Zq, after which it re- 
mains constant. The metal enrichment of stars in the 
fgas = 0.8 runs proceeds rapidly, increasing in a few 100 
Myr from 0.1 Zq to 0.5 Zq and leveling off after ~ 0.5 
Gyr at values between 0.8 Zq and Zq. 

The overall timespan covered by each simulation was 
2 Gyr. Fig. 1(a) illustrates the star formation history 
of one of the merger simulations. Specifically, we show 
a simulation starting with 40% gas fraction and ending 
with a stellar mass of 1.5 x 10^^ AIq, but the SFH of the 
other simulations looks qualitatively similar. Fig. 1(b) 
illustrates the build-up of stellar mass. In order to allow 
a fair test of the SED modeling with Bruzual & Char- 
lot (2003) templates (that include mass loss), we take 
into account mass loss of each of the stellar particles for 
which the time of formation and gas mass from which it 
was formed are stored in the simulation snapshots. This 
explains the slow decrease in total stellar mass at late 
times in Fig. 1(b). Fig. 1(c) shows the exhaustion of gas 
from which the stars are formed, and Fig. 1(d) presents 
the accretion history onto the black hole(s). We draw the 
time axis relative to the actual moment of merging, de- 
fined as the timestcp when the two black hole particles 
become one, coinciding with the peak in the accretion 
history. Cross symbols indicate the snapshots, separated 
by 70 Myr, when all physical information was stored to 
disk. 

As time progresses, the orderly rotation and star for- 
mation in the disks is disturbed by each others gravi- 
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Fig. 1. — Evolution of a typical merger simulation, (a) The 
star formation history, (b) the mass build-up, (c) the gas exhaus- 
tion, (d) the accretion rate history onto the black hole(s), (e) the 
evolution of the intrinsic (i.e., unattenuated) and attenuated V- 
band luminosity, (f) the distribution of effective visual extinctions 
(attenuated minus intrinsic V-band magnitude) corresponding to 
different viewing angles, and (g) the intrinsic and attenuated U — V 
color. Cross symbols in panel (e) and (g) present the evolution of 
the intrinsic photometry. For the attenuated photometry in panels 
(e) to (g), a binned representation is used where a darker inten- 
sity indicates a larger number of viewing angles. The dashed curve 
represents the median evolution. The dotted lines indicate the in- 
terval containing the central 68% of the viewing angles. The cross 
symbols mark the sampling of snapshots when the full physical 
information of all SPH particles was stored to disk. After a first 
bump in the star formation rate during the first passage of the 
progenitors, a peak in star formation is reached for a brief period 
during which several hundreds of solar masses of gas axe converted 
into stars. The typical extinction for a random line of sight is 
peaJiing around the same time. Shortly after, the accretion onto 
the supermassive black hole is maximal, coinciding with the merger 
between the two progenitor black holes. The reddest U-V colors 
are reached during the merger remnant phase. 
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tational pull. The star formation history shows a first, 
but rather shallow, bump during the first passage of the 
disks. Next, gravitational torques enable the gas to loose 
angular momentum and flow to the centers where it trig- 
gers a starburst (Larson & Tinsley 1978; Noguchi 1988; 
Hernquist 1989; Barnes & Hernquist 1991, 1996; Mi- 
hos & Hernquist 1994, 1996). Meanwhile, part of the 
inflowing gas is fed to the central supermassive black 
holes (SMBHs). Once the SMBHs grow massive enough, 
they produce a luminous quasar (Sanders et al. 1988a, b; 
Hernquist 1989; Sanders & Mirabel 1996; Genzel et 
al. 1998) whose feedback contributes to halting subse- 
quent star formation (Di Matteo et al. 2005; Springel 
et al 2005a), leaving a red spheroid galaxy as remnant 
(Robertson et al. 2006; Cox et al. 2006). 

2.2. Extracting photometry from the simulation output 

The evolutionary path as outlined in §2.1 is followed by 
the GADGET-2 code at a flne time resolution {At ~ 10^ 
yr). At sparser timesteps (70 Myr apart), the positions, 
masses, ages, and metallicities of all particles were stored. 
It is from these simulation snapshots that we derive the 
observed SEDs of the merger as a function of time. The 
applied technique is similar to that used by Robertson et 
al. (2007) to translate simulations of the most massive 
2 ~ 4 — 14 galaxies into observables. 

The light a virtual observer would receive from the sim- 
ulated merger, is composed of stellar and AGN emission, 
the latter only contributing significantly during a brief 
period of time. We ignore any contribution from emis- 
sion lines produced by the gas content of the galaxies, 
possibly contributing on the order of 0.1 mag in the op- 
tical broad-band photometry. Furthermore, we account 
for attenuation by interstellar dust and Lyman forest at- 
tenuation by the intervening medium between the red- 
shifted galaxy and the observer following Madau (1995). 
The combination of these steps, described in this section, 
leads to observables that are similar to the real obser- 
vations that we model with stellar population synthesis 
codes. 

First, we focus on the computation of intrinsic (i.e., 
unattenuated) magnitudes from the stellar component. 
Each of the stellar particles is treated as a single stellar 
population characterized by its mass, age, and metal- 
licity. We choose to use the Salpeter (1955) IMF, as 
was done in previous observational work (e.g. Forster 
Schreiber et al. 2004; Wuyts et al. 2007). We then in- 
terpolate the corresponding luminosity for each stellar 
particle from a grid of SSP templates with different ages 
and metallicities from the stellar population synthesis 
code by Bruzual & Chariot (2003, hereafter BC03). Fig. 
1(e) illustrates the evolution of the intrinsic rest-frame 
V-band luminosity for one of the simulations. 

For the AGN emission, we scale a luminosity- 
dependent template SED by the bolometric black hole 
luminosity given by the simulation. The template SED 
was derived from the optically blue (i.e., unreddened) 
quasar sample by Richards et al. (2006) with locally at- 
tenuated light being reprocessed as an IR bump longward 
of A > 1 iJ,m. A full discussion of the AGN template is 
presented by Hopkins, Richards, & Hernquist (2007). In 
most of our analysis, we will consider the stellar light 
only. §4.5 addresses the impact AGN can have on the 
outcome of SED modeling during the brief period when 



its contribution to the total light is significant. 

Galaxies, certainly in their actively star-forming 
phases, are not devoid of gas and dust. It is therefore 
crucial to account for the obscuring and reddening ef- 
fect dust has on the stellar and AGN emission. We 
compute the optical depth along the line of sight to- 
ward each stellar particle. To do so, we compiite the 
local gas density on a fine grid derived from the SPH 
formalism and the particle distribution (Hopkins et al. 
2005a) and integrate out from each particle along the 
line of sight to large distance. The simulations are based 
on the GADGET multi-phase ISM model developed by 
Springel & Hernquist (2003). This model calculates the 
local mass fraction in the warm/hot {T = 10" — 10^ K, 
diffuse, partially ionized) and cold (T = 10^ K, in the 
simulations both associated with molecular clouds and 
HI cloud cores) phases of dense gas, assuming pressure 
equilibrium between the two phases. Following Hopkins 
et al. (2005b), the attenuation along the line of sight 
is then derived from the density of the warm/hot-phase 
component only. Hopkins et al. (2005b) found typi- 
cally small volume filling factors (< 0.01) and cross sec- 
tions of the cold-phase "clumps" , motivating their ap- 
proach. The assumption that most of the lines of sight 
only pass through the warm/hot-phase component pro- 
vides effectively a lower limit on the optical depths. We 
use a gas-to-dust ratio equal to that of the Milky Way, 
{Ab/Niii)mw — 8.47 x 10^^^ cm^, with a linear scaling 
factor accounting for gas metallicities deviating from so- 
lar: Ab/Nui = (Z/0.02)(>1b/7Vhi)mw- As default, we 
adopt the Calzetti et al. (2000) attenuation law for the 
wavelength dependence of the optical depth. Changes 
in the synthetic photometry when adopting a SMC-like 
or Milky Way-like attenuation law from Pei (1992) will 
be discussed. The computation of optical depths was re- 
peated for 30 viewing angles, with directions uniformly 
spaced in solid angle d cos 9d(j). Fig. 1(f) presents the dis- 
tribution of effective visual extinction values (observed 
minus intrinsic V-band magnitude) as a function of time 
since the merger. The extinction varies in the follow- 
ing way. In the early stages typical extinction values 
are modest, with the exception for a few lines-of-sights 
were the disks are seen edge-on. The overall extinction 
along all lines-of-sight reaches a peak during the merger- 
triggered starburst and drops to very low values after 
star formation has ceased. 

Finally, in computing the observer-frame apparent 
magnitudes, we redshift the attenuated SED and con- 
volve it with the same set of filtercurves that we have 
deep observations for in the Chandra Deep Field South 
(CDFS; Wuyts et al. 2008): B435, Veoe, ^775, 2^850, J, 
H, Ks, [3.6/im], [4.5/im], [5.8/im], and [8.0/im]. Here, we 
apply the depression factors Da{z) and Db{z) given by 
Madau (1995) for the Lyman forest attenuation of the 
continuum between Lya and Ly/3 and between Ly/3 and 
the Lyman limit respectively. The flux blueward of the 
Lyman limit (Al = 912A) was set to 0, as is done by the 
HYPERZ code (vl.l, Bolzonella et al. 2000) that we use 
for SED modeling. 

In practice, it is computationally more convenient to 
interpolate the apparent magnitudes in a given passband 
for each of the stellar particles on a precompiled grid of 
BC03 apparent magnitudes at the redshift of interest. 
The internal dust attenuation is then applied using the 
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value of the Calzetti et al. (2000) attenuation law at the 
effective wavelength for that passband. We tested that 
this method, as opposed to attenuating the full resolution 
BC03 spectrum and then convolving with the filtercurve, 
leads to photometric differences of at most a few percent. 

We note that we never attempt to separate the light 
into the contribution from the two progenitors. Instead, 
we always study the total photometry, as if the merging 
system were unresolved. 

In future studies, we will employ full-fledged radiative 
transfer codes such as SUNRISE (Jonsson 2006) to trans- 
late the simulation snapshot information to observables. 
SUNRISE calculates scattering, absorption, and rcemis- 
sion of light passing through the cold and/or hot phase of 
the ISM, with the optional inclusion of a subgrid model 
that accounts for the obscuring effects of birth clouds in 
which young star-forming regions arc embedded (Groves 
et al. 2008). Recmission by dust only influences the spec- 
tral shape longward of the wavelength bands used in our 
SED modeUng. Preliminary analysis of SUNRISE SEDs 
indicates that the other effects might lead to a stronger 
extinction during the star-forming phases, while leaving 
the spheroid photometry unchanged with respect to the 
photometry computed in this paper. It will be interest- 
ing, in the light of the findings presented here, to investi- 
gate the shape of the effective attenuation law obtained 
with such a more sophisticated code, and its dependence 
on the parameters of the radiative transfer. For the sake 
of this paper, our line-of-sight photometry has the ad- 
vantage of being transparent and simple. 

2.3. The color evolution of merging galaxies 

With the synthetic photometry at hand, we are now 
able to follow the color evolution of simulated galaxies 
throughout the merger event. We illustrate the spa- 
tially resolved color evolution for a face-on view of a 
typical simulation in Fig. 2. This particular simulation 
has an initial gas fraction of 80% and final stellar mass 
of 1.2 X 10"'^^ Mq. The 3-color images are composed of 
rest-frame U (blue), V {green), and J (red). On each 
postage stamp, we indicate the snapshot number, time 
since the merger, and whether the intrinsic (Int) or at- 
tenuated (Att) colors are plotted. The interval between 
consecutive snapshots is 70 Myr. The color distribution 
is clearly not homogeneous throughout the merger evo- 
lution. Instead, color gradients are present, with the 
nuclear regions being intrinsically bluer due to new star 
formation, but effectively redder than the surrounding 
material due to the presence of dust. We return to this 
point in an analysis of observed versus simulated galaxy 
colors by Wuyts et al. (2009b). In this paper, we restrict 
ourselves to the study of integrated colors. 

The corresponding evolutionary tracks in integrated 
rest-frame U — V versus V — J color-color space are pre- 
sented in Fig. 3. Labbe et al. (2005) first introduced the 
observed-frame equivalent of this diagram to illustrate 
the wide range of galaxy types at high redshift ranging 
from blue, relatively imobscured star-forming systems to 
dusty starbursts to quiescent red galaxies. The color cri- 
terion proposed by Labbe et al. (in preparation) to select 
galaxies of quiescent nature is plotted as the red wedge 
in Fig. 3. 

Considering the intrinsic (unattenuated) colors, the 
evolutionary track has a direct relation to the mean light- 



weighted age (and depends to a lesser degree on the grad- 
ual increase in stellar metallicities). It bends down to 
bluer U—V as new star formation is triggered during first 
passage of the progenitors, and again during the shorter- 
lived nuclear starburst phase. Eventually, the color track 
leads to the region in color-color space where also quies- 
cent observed galaxies reside. In reality, the path taken 
in color-color space will not merely depend on the star 
formation history. Taking into account attenuation by 
dust, the integrated colors get redder. This is particu- 
larly true during the phases of active star formation. We 
note that the precise shape of the color track depends 
not only on initial conditions (gas fraction, mass, disk 
orientations, ...) but also on viewing angle. In fact, as 
illustrated in Fig. 3, color differences between different 
viewing angles at a given time can be of similar size as 
color changes between different times for a given viewing 
angle. At later times, the role of dust is minimal and the 
remnant converges to colors typical for quiescent galaxies 
(i.e., within the red wedge). 

2.4. The colors and SEDs of simulated and observed 
galaxies 

Prior to analyzing the performance of our SED mod- 
eling procedure, it is important to confirm that the sim- 
ulated galaxies have spectral shapes resembling those of 
real high-rcdshift galaxies in observed deep fields, thus 
validating their role as test objects. To this end, we indi- 
cate the binned color distribution of simulated galaxies, 
viewed from different angles and during different phases 
of their evolution, in a rest-frame U — V versus V — J 
color-color diagram (Fig. 4). Plus symbols show the lo- 
cation of observed galaxies in the HDFS (Labbe et al. 
2003), MS 1054-03 (Forster Schreiber et al. 2006), and 
the CDFS (Wuyts et al. 2008) selected by their photo- 
metric redshift (or spectroscopic when available) to lie in 
the same redshift range (1.5 < z < 3.0). We also applied 
a stellar mass cut at M^, > 1.4 x 10^" Mq for the observed 
sample; the lowest initial stellar mass for the considered 
set of simulations. Here, we do not attempt to statisti- 
cally compare the two samples. The abundances of dif- 
ferent types of galaxies as predicted from the simulations 
will be addressed by Wuyts et al. (2009b). For our cur- 
rent purpose of analyzing the effects of star formation his- 
tory, dust, metallicity and AGN on SED modeling, it is 
sufficient to note that there is a large overlap between the 
color-color space spanned by the simulated and observed 
galaxies. Every simulated galaxy considered in this pa- 
per has counterparts in the real universe with similar 
rest-frame optical and optical-to-NIR colors. However, 
the observed distribution extends to colors that are red- 
der by a few 0.1 mag, both in U—V and in V—J, than the 
simulations considered here. Given the one-sided nature 
of the different color spread, it is unlikely that this can 
be attributed to photometric uncertainties alone. There- 
fore, we caution that our results may not necessarily be 
extrapolated to the reddest galaxies present in observed 
samples. Wuyts et al. (2009b) discuss several possible 
origins of the color discrepanc;y, particularly in V — J, 
ranging from the method to compute model photometry 
to differences between the evolutionary history of high- 
redshift galaxies in the real universe and in our merger 
simulations. There we analyze, for example, the pres- 
ence of color gradients, and the dependence of the color 
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Fig. 2. — Three-color postage stamps (rest-frame U (blue), V (green), J (red)) of a typical merger simulation. The top three rows 
illustrate the intrinsic evolution, ignoring extinction. The bottom three rows show the evolution when taking into account attenuation by 
dust. Significant changes in color are observed over time, and between intrinsic and attenuated images. 



distribution on the adopted attenuation law and stel- 
lar population synthesis code. Other models than BC03 
that include a more significant contribution of the ther- 
mally pulsing asymptotic giant branch (particularly in 
the near-IR) result in a color distribution that is shifted 
towards redder V — J colors (Maraston 2005; Maraston 



et al. 2006). 

To ascertain that observed and simulated galaxies with 
similar U—V and V~J colors have similar SEDs over the 
whole spectral range, Fig. 5 presents the rest-frame SEDs 
of objects in region 1-6 of Fig. 4. Again, the binned distri- 
bution represents the simulations, with darker grayscales 
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Fig. 3. — Evolutionary tracks in a rest-frame U — V versus V — J color-color diagram for the typical merger simulation shown in Fig. 
2. We illustrate the intrinsic color evolution as well as the evolution of attenuated colors with time as we view the merging pair face- and 
edge-on. Snapshots 70 Myr apart are marked with small dots. Key phases corresponding to the numbers indicated in Fig. 2 and Fig. 14 
are indicated with large symbols. Time since (or before) then actual moment of merging is indicated in the top right corner of each panel. 
Colors can vary significantly at a given time (by ~ 0.5 mag) as we view the merging system from different lines of sight. The starburst 
occuring during the merger is very dusty, leading to redder colors at that stage. 



2.5 




-0.5 r 1 

-0.5 0.0 0.5 1.0 1.5 2,0 2.5 

(V-J).est 

Fig. 4. — Rest-frame U — V versus V — J color-color diagram 
showing the binned color distribution of the simulations seen under 
different viewing angles and at different epochs. Overplotted {plus 
symbols) are the rest-frame colors of observed galaxies with A/* > 
1.4 X 10^'' Mq at 1.5 < z < 3 in the HDFS, MS 1054-03, and the 
CDFS. Observed galaxies with matching colors are found for all 
simulated galaxies. The reddest observed sources in U — V and 
V — J are not reproduced by the considered set of simulations. 
Rest-frame SEDs for sources in regions 1-6 are displayed in Fig. 5. 



indicating a larger number of objects. Overplotted with 
black dots is the broad-band photometry of our ob- 
served sample within the same region of color-color space, 
placed at the respective rest-frame wavelength. The 
SEDs are normalized to the rest-frame F-band. By se- 
lection, the observed and simulated photometry matches 
well at rest- frame U and J. In between the UVJ filters, 
and outside the U-to-J range, no correspondence was im- 
posed. The fact that the UV spectral shape and the NIR 
tail of the observed and simulated SEDs show a general 
agreement, is encouraging. We conclude that the simu- 
lated photometry can be adopted as a realistic input to 
our SED modeling procedure. The results of our anal- 
ysis will be applicable to observed galaxies with similar 
colors. 

3. SED MODELING: METHODOLOGY 

We characterize physical parameters such as stellar 
mass, stellar age, dust attenuation, and SFR by matching 
the observed-frame broad-band photometry to synthetic 
templates from the stellar population synthesis code by 
BC03. We use the HYPERZ stellar population fitting 
code, version 1.1 (Bolzonella et al. 2000) and fit the SED 
twice: first fixing the rcdshift to the true value (for which 
we computed the simulated photometry), next adopt- 
ing a photometric redshift estimate obtained from the 
EAZY version 1.0 photometric redshift code (Brammer 
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Fig. 5. — Rest-frame SEDs of simulated galaxies in regions 1-6 of Fig. 4. A darker intensity of the binned representation indicates a larger 
density of simulated galaxies with that flux level. In each panel, we give the central 68% interval of the distribution of times before or since 
the merger for the simulation snapshots with photometry in the respective region. Overplotted {black dots) are the rest-frame broad-band 
SEDs of observed 1.5 < z < 3 galaxies with M, > 1.4 x 10^° Mq in the HDFS, MS 1054-03, and the CDFS. A general agreement between 
observed and simulated spectral shapes is observed, also outside the U-to-J range where the correspondence was not imposed by selection. 



et al. in preparation). In each case, the full B435-to-8 
/xm SED, sampled with identical passbands as available 
for the GOODS-CDFS (B435, T^eoe, «775, ^sso, J, H, K.,, 
[3.6 y^m], [4.5 /im], [5.8 /iin], [8.0 /xm]), was fed to HY- 
PERZ. Random photometric uncertainties were assigned 
as to mimic real observations in the CDFS, and fluxes 
in each band were perturbed accordingly. Precisely, for 
each of the 43200 SEDs corresponding to a simidatcd 
galaxy observed during a certain phase of its evolution, 
placed at a certain redshift, and observed along a cer- 
tain line-of-sight, wc compute 5 realizations of the SED 
by introducing a gaussian perturbation in all bands with 
tiic ampUtudc derived from the depth of GOODS-CDFS 
observations in the respective bands. A minimum error 
of 0.08 mag was adopted for all bands, preventing small 
errors from dominating the fit. 

As in Wuyts et al. (2007), we selected the least 
solution out of three possible star formation histories: a 
single stellar population (SSP) without dust, a constant 
star formation (CSF) history with dust {Ay varying from 
to 4 in steps of 0.2), and an exponentially declining 
star formation history with a fixed e-folding timescale of 
300 Myr (T300) and identical range of Ay values. Ages 
were constrained to be larger than 50 Myr, to prevent 
improbably young ages, and smaller than the age of the 
universe at the observed redshift. We used a Calzetti 
et al. (2000) attenuation law, and assumed solar metal- 
licity and a Salpcter (1955) IMF with lower and upper 
mass cut-offs Q.IMq and lOOM©. In recent literature. 



several authors have used an IMF with fewer low mass 
stars, such as presented by Kroupa (2001). We note that 
the choice of IMF is not the focus of this study, and our 
results remain valid as long as we use the same IMF to 
compute the synthetic photometry as to model the re- 
sulting SEDs. Likewise, we use for consistency the same 
stellar population synthesis code (BC03) to compute and 
to fit the synthetic SEDs. 

When refering to the age derived from SED modeling, 
we mean the mass-weighted age obtained by integrating 
over the different ages of SSPs that build up the best-fit 
SFH, weighted with their mass fraction taking into ac- 
count mass loss over time. This measure aims to quantify 
the age of the bulk of the stars. For an SSP, it equals the 
time passed since the single burst. For a CSF history, 
it is essentially half the time passed since the onset of 
star formation. The T300 SFH represents an intermedi- 
ate case. 

4. RESULTS FROM SED MODELING AT FIXED REDSHIFT 

First, in §4.1, we present the overall performance of 
the standard SED modeling applied to the 'full' photom- 
etry, taking into account the effects of both attenuation, 
mctallicity, and AGN contribution as realistically as pos- 
sibles Sec;ond, we discuss the impact from different as- 
pects influencing the colors and luminosities of galaxies 
separately. In order to isolate effects from star formation 
history (§4.2), dust attenuation (§4.3), metallicity varia- 
tions (§4.4), and AGN contribution (§4.5), we computed 
the photometry for each snapshot with and without at- 
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tenuation, with and without AGN contribution, and us- 
ing solar mctaUicity, or the mctaUicity as computed by 
the simulation for each stellar particle. To each of these 
sets of SEDs, we applied the modeling described in §3. 
This approach enables us to reconstruct the analysis step 
by step, adding one aspect at a time. Finally, we discuss 
what our analysis means for galaxies with different rest- 
frame colors (§4.6). 

4.1. Overall performance 

The combined effects of mismatch between true and 
template SFH, attenuation by dust, mctaUicity varia- 
tions and AGN activity on our ability to characterize the 
stellar population of a galaxy are summarized in Fig. 6. 
Here, we plot the best estimate from SED modeling ver- 
sus the true^ value of the considered stellar population 
property. Wc bin the distribution of points for different 
initial conditions, timesteps and lines-of-sight. Darker 
intensities represent a higher density in the bin. The 
solid curve represents the median of the distribution and 
the dotted curves mark the central 68% interval. 

The youngest ages in Fig. 6(a) correspond to the early 
stages of simulations with an initial gas fraction of 80%. 
The stars present at the start of these simulations were 
assigned young ages and low metallicities (see §2.1). In 
§4.4, we describe how modeling sub-solar metallicity pop- 
ulations with solar metallicity templates can lead to age 
overestimates. Since the lower gas fraction runs start 
with higher stellar ages and metallicities, they do not 
show this trend. Fig. 6(b) compares the recovered and 
true stellar mass, normalized to the final stellar mass 
of the simulation. At low M/Mfinai ratios, i.e., at the 
start of the high gas fraction simulations, the mass esti- 
mates agree well with the true values. H(^rc;. the rciason 
is that biases in age and extinction average out. On the 
one hand, the overestimated age leads to an overesti- 
mate of the intrinsic mass-to-light ratio ( 



M 



Lin 



> 



M 



). On the other hand, the amount of extinction 



tends to be underestimated ( 



< 



LAt 



)• 



As a result, the mass estimate for this phase is relatively 
robust: 



estimate 



= LAt 



L 



Int 



^Att 



Lint 



M 



estimate 

M 



L 



Int 



(4) 



^Int 



(5) 

true 

(6) 

is the attenuated 



= Mtrue 

where Lint is the intrinsic and ^Att 
(i.e., observed) luminosity. The largest systematic mass 
underestimates occur at intermediate M/Mfi„ai, during 
the merger-triggered star-forming phases of tlw^ simula- 
tion. This can be explained in terms of a mismatch be- 
tween the true SFH and the SFH of the best-fitting tem- 
plate (see §4.2), in combination with an underestimate of 

^ Hereafter, we use 'true' when referring to the value of a stellar 
population property as computed using information on all SPH 
particles in the simulation. Quantities derived by fitting templates 
to integrated SEDs are referred to as 'estimates'. While 'true' mass 
and mass-weighted age can be read directly from the simulation, 
wc note that parameters as 'true' reddening and visual extinction 
depend on the adopted stellar population synthesis code and details 
of the radiative transfer, which themselves may have biases. 



the extinction (see §4.3). Briefly, such a SFH mismatch 
leads to age underestimates, hence underestimates, 
and, added to an insufficient dust correction, to mass 
underestimates. Finally, the correspondence is best at 
log(M/Mfinai) ^ 0, where the merger remnants reside. 

The reddening (Fig. 6(c)) is overall well reproduced. 
Only at the highest reddening levels, the agreement de- 
teriorates. The latter correspond to the times when and 
viewing angles under which the effect of increased ex- 
tinction toward young stars is maximal, as will be ex- 
plained in §4.3. As opposed to the reddening, however, 
the extinction (Fig. 6(d)) shows large systematic under- 
estimates. Using 

Air 

Ry = = 4.05 (7) 



E{B - V) 



to translate the selective absorption E{B — V) into a 

total visual absorption Ay results in an underestimate 
over the whole range of Ay values, particularly during 
the highly obscured phases. An observer is limited by the 
light that he/she receives. In §4.3, we discuss in depth 
how the sum of emitting sources that are each attenuated 
according to Calzetti et al. (2000) does not follow that 
same reddening law. 

Finally, since the reliability of SFR indicators from X- 
rays over optical and IR to radio wavelengths has been 
topic of much recent debate (e.g., Reddy et al. 2006; 
Papovich et al. 2006; 2007; Daddi et al. 2007a), it is 
interesting to investigate its recovery by modeling the 
optical-to-8 /im SEDs. This is presented in Fig. 6(e). At 
low and intermediate SFRs {SFRtmc < 100 Mq/jt), we 
find a satisfying agreement, with a median offset that in- 
creases slightly towards higher SFRs. For more than 84% 
of the galaxies with SFRnuc > 100 Mq/jt, we underes- 
timate the SFR, typically by factor of 3. These systems 
often suffer strong dust obscuration, and therefore this 
bias in SFR directly relates to the large underestimates 
of extinction at the high Ay end in Fig. 6(d). 

In order to quantify the performance of the SED mod- 
eling in estimating stellar mass-weighted ages, wc define 
Alog{agew) as log{ageni,recovered) - logiageu,,true)- Sim- 
ilar definitions are used to quantify the offset in mass, 
reddening and extinction, always indicating an underes- 
timate with a negative value of A = parameter recovered- 
parameter true- Fig. 7 presents the performance of the 
SED modeling on the full photometry (including dust, 
metallicity variations, and AGN), expressed by the A 
values, as a function of time since the merger. We adopt 
the same binned plot style as in Fig. 6, with darker inten- 
sities meaning a higher density in the bin. Open boxes 
contain less than 1% of the total number of SEDs at that 
timestep. 

We quantify the performance of the SED modeling sep- 
arately for galaxies in the 'disk', 'merger', and 'spheroid' 

regime by computing the median and the central 68% 
interval of the distribution of A values for all sim- 
ulation snapshots (under a range of viewing angles) 



in that phase. We find Alogage^ 



,disk 



-h0.03 



+0.19 
-0.42' 



Alogage^^merger = -0-12Iq;26, and Alogage^^spheroid = 
—0.03^0^4. The underestimate and scatter is largest for 
the phases of merger-triggered star formation. These are 
the statistics for low and high gas fraction runs com- 
bined. Considering simulations with an initial gas frac- 



10 



o 



CD 
O 

Ol 
O 



0.5 



0.0 



-0.5 



1.0 



■1.5 



III, 


111,111,1 

□ □ □ .□ 

Lj _ _ n 


1 , 1 1 1 1 1 1 1 , 1 1 1 , 1 

□ □ □ 1 

□ u. u MMm ■ 

■ 

^•SM* ■. ■■ 




. . ■ ■ ■ ■ 
■ p'i ■ ■ □ 






' y'' ■ ■ ^ ■ 


□ ■ 


; □□□□□□□ 
J □ □ □ □ □ □ 
□ □ 

\ ID □□ 
■■■ \ 1 □ □ □ □ 


□ □□□□ 

□ □ □□□□ 
□ □□□□□□ 

□□ □□□□□□□□ 
□□□□□□□□□□□ 
□□□□ □□□□□□ 
□ □ □□□□□ 




... 1 


. . . 1 . . . 1 . 


. 1 ... 1 ... 1 ... 1 . 



o 

E 



CO 

m 

D 

o 



0.5 



0.0 



-0.5 



1.0 



■1.5 



1 

■ (b) 


1 






1 










-1 r-i H 










^.--^ / ■ 










. .-■ f' 










'i£il&IIPrJig;Hluriri 




s 'a' 










!!! «■ 






■i*L i?.(i_'!Ji J 




'oi 'i^ if 










^. ! < ^ 






'SI:W)«ti . 1 i"i i II I 




111 «a \Zl 




M i 




ssr 








£.._'!. II.J 

^■"ii'iri 


Hi 

1 


□ I, 

1 




J 


1 



-1.0-0.8-0.6-0.4-0.2 0.0 0.2 0.4 
log(age„_i,J [Gyr] 



-0.6 -0.4 -0.2 0.0 
log(masSi,„^ / nnass,r„„,,i,J 




0.0 



0.2 0.4 
(B-VU - 



0.6 
(B-V)„i 
4 



0.8 



. 1 1 ' ' ' 

; (d) 


1 1 1 l_. 1 , 1 1 1., 1 


' ' ' ' 1 y 


□ □□ □ 






□ □□□ □ 


□ 




□□□□ □ □ 






□□□□□□□□ 


□ 






n - - / 




I": u li l: :; □ v.- :, .: l 






' J J i_ _: □ n !□ ■ i: 


mum / 




G c □ ■ ■ ■ ■ ■ 


■ ■ mmti _ _ 




1 1 M ■■P"«B 


■ ■^^■■•■-■■■■B 




□ □ i: ■ r 






11 1 n B 1 / 




■ a et u V! 


/f ■ 










J ■ ■ ■ Hta 






■ III \~A 














f ^'K', ^ ¥ ^ F.', F 




V . . . 1 " 


1 


2 


3 ^ 



^ 3 

© 



Li_ 
O 



.1 

; 


' 1 ' 


1 1 1 








■ / 
















/-■■■»■ 








/ ■ 
















JGOD X*'^""^ 
























' □ a a a a'a a in >n a a'a a a 
















'i a a a.a aauaaaa a c □ a ^ 








'aMaaM">*>" ''"j 








..■^ a 








fiaX X"^'~ 








.■ / ipwn 








/ X'J"'' •■ 








aX^Aaaaa J iaanu 








/ AuwM a 








^munu u □□□□ f- 
















1 .laaBGDaDDDDD l: 








/■ BBBnggnn □ d 






■ .-■ 






- ■■/ 


»■■ 






: ' 


■ 








Fi-i^ 


. . . 1- 1 





-10 1 2 3 4 

log(SFR,,J [Mo/yr] 

Fig. 6. — Overall performance of the SED modeling. Recovered versus true (a) mass-weighted stellcir age, (b) ratio of current to final 
stellar mass, (c) effective reddening (i.e., attenuated minus intrinsic B — V color), (d) effective visual extinction (i.e., attenuated minus 
intrinsic V-band magnitude), and (e) star formation rate. The SED modeling was performed on the total (stellar+AGN) attenuated 
photometry corresponding to simulations with a range of masses, starting with 40% and 80% gas fractions, and seen at different timesteps 
and under different viewing angles. The solid curves indicate the median and dotted curves comprise the central 68% of the distribution. 
The total visual extinction Ay is the least constrained of the studied parameters. In particular for heavily extincted galaxies the Ay is 
greatly underestimated. 
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Fig. 7. — Overall performance of the SED modeling. The difference between estimated and true (a) mass-weighted age, (b) stellar mass, 
(c) effective reddening, (d) effective visual extinction, and (e) star formation rate as a function of time since the merger. The SED modeling 
was performed on the total (stcUar+AGN) attenuated photometry. The solid curves indicate the median and dotted lines comprise the 
central 68% of the distribution. The properties of merger remnants are well rcijroduccd. The results for star-forming galaxies, especially 
for those in the phase of merging, show underestimates in both age, extinction, star formation rate, and mass. 
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tion of 40% separately, the age underestimates during 
the star- forming phases increase: Alogage^^disk,4o% = 
-0.13+°:^^ and Aloga5e^,^erger,40% = -O-SOtHi- The 
/gas = 80% runs instead show age overestimates until 
^ 200 Myr before the merger, mainly due to their lower 
initial metallicity (see §4.4). During the spheroid phase, 
ages are recovered similarly well for the low and high 
gas fraction runs. The recovery of other stellar popula- 
tion properties does not depend strongly on gas fraction 
unless stated otherwise. 

Fig. 7(b) shows the same scenario as described for Fig. 
6(b). Namely, mass estimates are robust once the merger 
remnant phase is reached (AlogMsphcroid = ^0.02tQ;'j'i), 
but more significant mass underestimates occur for star- 
forming disks (AlogAfdisk — ^0.06^q''|'4); and espe- 
cially during phases of merger-enhanced star formation 

(AlogM^erger = -0.13t°;l°). 

Averaged over all SEDs for a given time with re- 
spect to the merger (representing different viewing an- 
gles, masses, and gas fractions), our reddening esti- 
mates are well behaved: AE{B - V')disk = -0.02tam, 
A^(B-y),„ergcr = -0.02t|];°«, and AE{B-V)spbe.oid = 
—O.OSf^l}. However, as noted already in Fig. 6(d), 
the bias in estimated visual extinction is often much 
larger, particularly during the star-forming (disk and 
merger) phases: AAy,disk = -0. 35^0.34, AAy^merger = 
-0.54t°;^°, and AAy,3phcroid = -0.29t° j^. Part of the 
scatter quoted for AE{B — V) and AAv results from dif- 
ferences bctwecin the low and high gas fraction riuis. The 
three key differences, with reference to where their im- 
pact is detailed, are the following. First, the high gas 
fraction runs start out with lower metallicities (§4.4). 
Second, the low gas fraction runs start with older stars. 
The larger age range makes the effect of age-dependent 
extinction more pronounced (§4.3.2). Finally, the high 
gas fraction runs feed relatively more gas to the SMBH, 
resulting in a larger AGN contribution to the integrated 
light during final coalescence (§4.5). Since the obscura- 
tion during the merger phase tends to be to the newly 
formed stars mainly, the underestimated extinction leads 
to an underestimate of the SFR. Consequently, the evo- 
lution of Alog SFR with time is as a mirrored version 
of the star formation history, showing the largest un- 
derestimates during first passage and final coalescence. 
We find Alog 5Fi?disk = -0.22t°;i, A log 5i^E^erger = 
-0.44t°:32, and A log ^^^i^sphoroid = -0.231°;^?. Note 
that the offset for the spheroid phase merely means that 
we arc unable to distinguish low (a few 0.1 MQ/yr) from 
very low ('^O.l Mq/vy) star- forming systems. We are 
able to robustly classify them as quiescent. 

In the above statistics, each A value corresponding to 
a particular snapshot was given equal weight. Alter- 
natively, we can quantify the performance of our SED 
modeling in a global sense by summing the SFRs of 
all snapshots and comparing this to the sum of all re- 
covered SFRs. Such an approach gives more weight to 
phases with higher SFR, as is also the case for observa- 
tional studies of the cosmic SFR density. Considering 
the entire timespan for all simulations, we find that our 
SED modeling underpredicts the total mass of stars being 
formed per unit time by a factor 1.8. A similar exercise 
shows that the total assembled stellar mass is underes- 
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Fig. 8. — Correlations between biases in the estimated age, mass, 
reddening, extinction, and SFR obtained from broad-band SED 
modeling. Style as in Fig. 7. Positive values of A mean overes- 
timates, negative values refer to underestimates with respect to 
the true age, mass, E{B — V)ea, Ay^eS, and SFR in the simu- 
lations. Strong correlations exist between A values for different 
stellar population properties. E.g., over/underestimates in redden- 
ing E{B — V)gff lead to offsets of the same sign in extinction Ay eff 
and SFR. 



timated by only 14%. When focussing on the merger 
regime alone, the underestimates increase to a factor 2 
for the summed SFR and 30% for the summed stellar 

mass. 

We present an overview of correlations between the A 
values of all considered properties in Fig. 8. It is clear 
that biases in estimates of stellar population properties 
are often strongly correlated. For example, when the 
age is severely underestimated, the same will be true for 
stellar mass. Likewise, a strong positive correlation ex- 
ists between biases in the estimated reddening, extinc- 
tion, and star formation rate. In the following sections, 
we shed more light on how the estimates of age, mass, 
reddening, extinction, and SFR are connected. 

In the light of the upcoming Spitzer warm mission, 
we note that the interval containing the central 68% of 
AlogM values becomes twice as broad when omitting 
the IRAC bands. This test demonstrates the value of 
extending SEDs to longer wavelengths. Similar improve- 
ments by including IRAC photometry are found for the 
other stellar population parameters considered in this pa- 
per. 

We performed our analysis on synthetic photometry of 
galaxies placed at redshifts z = 1.5 to z = 2.9. We note 
that the results do not change if you change the redshift 
of the merger. Furthermore, the trends described in this 
section are all systematic and cannot be attributed to 
signal-to-noise effects (e.g. more cxtincted galaxies at 
the highest redshifts being fainter and therefore less well 
recovered). We tested this in two ways: first by omitting 
the perturbation of the synthetic fluxes, second by ap- 
plying a conservative cut in the observed iiTg-band mag- 
nitude: Ks^ohs < 23.6, corresponding to S/Nk^ > 10. In 
both cases, the same trends described in this Section are 
still present. 
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4.2. Impact of mismatch between true and template SFH 

The biases in estimates of stellar population proper- 
ties (§4.1) result from a complex interplay between SFH, 
dust distribution, metallicity variations, and AGN con- 
tribution. Wc now separate these effects in order to ob- 
tain a better understanding of the origin of the biases. In 
this section, we repeat the SED modeling on the intrin- 
sic stellar photometry (i.e., no dust extinction or AGN 
contribution), with the metallicity of the stars fixed to 
solar (i.e., no metallicity variations). Each stellar par- 
ticle is still treated as a SSP with an age given by the 
GADGET-2 output (i.e., the SFH is identical to that on 
which the 'full' photometry in §4.1 was based), but we 
extract the photometry from solar metallicity SSP mod- 
els, even if the stellar particles were recorded to have 
lower/higher mctallicities in the GADGET-2 simulation. 
Since we use solar metallicity templates to compute the 
integrated light of the simulations as well as to fit the 
resulting synthetic SEDs, this procedure removes effects 
due to mismatch in metallicity and isolates the role of 
the SFH. This illustrates the consequences of a possi- 
ble mismatch between the true SFH in the simulation 
and the allowed template SFHs in our SED modeling. 
Our aim is to constrain the properties of the bulk of the 
stars (e.g., mass- weighted age), whereas our input SEDs 
are obviously of a light-weighted nature. Since massive 
O and B stars make young stellar populations brighter 
than older stellar populations, giving them more weight 
in the integrated SED, the light-weighted stellar age will 
be younger than the mass-weighted stellar age. This is 
always the case, but provided we have a template rep- 
resenting the correct SFH, it is possible to account for 
this effect and still find the correct age of the bulk of 
the stars. Our three allowed SFHs are an SSP, where all 
stars formed in a single burst, a model with SFR oc e"*/"^ 
with T = 300 Myr, and a constant star formation history. 
These are standard choices in analyses of distant galax- 
ies. However, they do not encompass a star formation 
history where the rate of star formation was lower in 
the past than it is now, as is the case during first pas- 
sage and during the actual merger-triggered starburst 
(see Fig. 1). In general, fitting a template SFH that has 
< l^^^^ltrupi older population will 
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be lost to some degree under the glare of newly formed 
stars, leading to an underestimate of the age. 

This is illustrated in Fig. 9(a), where the difference 
between recovered and true mass-weighted stellar age is 
plotted as a function of time since the merger between 
the SMBHs. During the first snapshot, when the star for- 
mation history matches (by construction) our CSF tem- 
plate, we find no systematic offset and a low scatter, 
purely resulting from photometric uncertainties. Soon 
after, we start to underestimate the age, with minima 
coinciding with the moment of first passage (500 Myr 
before the actual merger) and that of the actual merger- 
triggered starburst. It is precisely at these moments that 
the real SFH deviates most from the allowed template 
SFHs (during the merger, the t^qq and CSF templates are 
most frequently preferred by the fitting routine) . During 
the starburst phase itself, the median offset of true mass- 
weighted age versus recovered age exceeds 0.6 dex, with 
a large scatter due to differences in the SFH for different 
initial conditions. For example, the ratio of SFR at first 



1.0 - (o) 



0.5 





en 
D 

cn 
O 




■1.0 



■1.5 



0.5 



0.0 



m 
If) 
o 

o 



-0.5 



■1.0 



■1.5 



Li_ 
CO 



o 



ldisk merger spheroid 



-0.5 0.0 0.5 1.0 1.5 

Time since merger [Gyr] 



(b) 



□ □ □ n 



.DISK MERGER SPHEROID 



-0.5 0.0 0.5 1.0 1.5 
Time since merger [Gyr] 



(c) 



m m /'mmm 

"'"^ 0. vm 




DISK MERGER 



SPHEROID 



-0.5 0.0 0.5 1.0 

Time since merger [Gyr] 



1.5 



Fk;. 9. — Impact of star formation history. The difference bc- 
t"ween estimated and true (a) mass-weighted age, (b) mass, and (c) 
SFR as a function of time for all simulations, with the SED mod- 
eling performed on the intrinsic (i.e., unatteiiuated) stellar pho- 
tometry -with all stars set to solar metallicity. Style as in Fig. 7. 
Deviations from (negative indicating an underestimate) are due 
to mismatch between the actual star formation history and the his- 
tories allowed in our SED modeling (SSP/CSF-|-dust/T3oo+dust). 
Maximal underestimates of age and mass are reached during the 
merger itself. A secondary minimum is reached during first passage 
of the progenitors, ~ 0.4 Gyr before. The SFR, however, is well 
recovered during the episodes of merger-enhanced star formation, 
but not in the first few 100 Myr after. 
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passage over SFR during the central starburst increases 
with gas fraction. After aU activity has quieted down, 
the derived ages he within 0.1 dex of their true value. 

Since one tends to count the young light only, mass 
will be underestimated as well during first passage and 
the final nuclear starburst (see Fig. 9(b)). For the same 
reason, models allowing for a secondary burst of star for- 
mation on top of an older stellar population were found 
to reveal larger total stellar masses, in particular for blue 
objects (Papovich et al. 2006; Erb et al. 2006b; Wuyts 
et al. 2007). As for the age measurement, the derived 
masses for the merger remnants lie within 0.1 dex of their 
true value. 

Since no knowledge about underlying old stellar popu- 
lations is needed to constrain the ongoing star formation 
rate, SED modeling of the intrinsic light reproduces the 
SFR well during the episodes of enhanced star formation 
(Fig. 9(c)), when age and mass determinations arc least 
reliable. However, the effect of template mismatch be- 
tween true and best-fit SFH comes into play later. In 
the simulations, the SFR drops rapidly after reaching its 
peak during final coalescence. The preferred SFH of the 
best-fit model most frequently has a slower decrease, over 
an e- folding time of 300 Myr. By accounting for the bulk 
of recent star formation, the model consequently overpre- 
dicts the ongoing SFR during the first few 100 Myr after 
the merger, by ^ 0.6 dex. Once more time has passed, 
the SED modeling correctly identifies the remnant as a 
low SFR galaxy. 

We check that the systematic offsets and scatter in Fig. 
9 are indeed due to mismatches between the true and 
template SFHs by performing the following test. For 
each simulation, we change the distribution of ages of 
stellar particles to reflect the assimied SFHs that are 
used in the SED fitting (see §3). We then compute the 
synthetic photometry in the redshifted passbands as de- 
scribed in §2.2, still ignoring dust reddening and assum- 
ing solar metallicity. Modeling the resulting idealized 
synthetic SEDs, we recover the mass- weighted age, mass, 
and SFR essentially without systematic offsets (< 1%), 
and with little scatter (a fifth of the scatter obtained 
when the true SFH from the simulations was adopted). 
The scatter in A log M decreases by only 40% when omit- 
ting the perturbation of the synthetic fluxes, and for 
Alogageto and IS.\ogSFR the change is even smaller. 
This indicates that the discrete sampling of the SED by 
the observed passbands and the discreteness of the age 
steps allowed in the flt contribute by a similar amount 
as photometric uncertainties to the very small A values. 
This test demonstrates that, for a survey as GOODS, 
mismatch between true and template SFH is a much 
larger limitation to constrain stellar population prop- 
erties of 1.5 < z < 3 galaxies than, e.g., photometric 
uncertainties. 

4.3. Impact of attenuation 

4.3.1. Non-uniform attenuation uncorrelated to 

properties of emitters 

As illustrated in §2, dust attenuation has a large effect 
on colors, SEDs, and fluxes in all star- forming phases. 
This is the case for both isolated spiral galaxies (see also 
Rocha & Jonsson 2008) and interacting systems ( Jonsson 
et al. 2006). As described in §3, we use the approach of 
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Fig. 10. — The effective reddening (attenuated minus intrin- 
sic B — V coior) versus total absorption in the V-band for all 
tinicstcps, viewing angles and initial conditions. The intensity of 
the bimied distribution indicates the number of simulations in the 
respective part of the diagram. A ratio of total to selective ab- 
sorption Rv = 4.05 as by Calzctti et al. (2000) is plotted with 
the thick dashed line. The dot-dashed curve indicates a series of 
toy models where the distribution of Ay values is uniform between 
and a maximum value, and all emitting sources are identical. 
Stellar particles individually have Rv = 4.05, but in the case of a 
non-uniform dust distribution the sum of all stellar particles has 
an effective Rv > 4.05. 



a uniform foreground screen to account for the attenua- 
tion by dust in our SED modeling. Fig. 1(f), illustrating 
the range of effective visual extinction values (attenu- 
ated minus intrinsic y-band magnitude) for a random 
simulation depending on the viewing angle, proves that 
such a representation is not valid. Here we address the 
impact that a non-uniform distribution of the dust will 
have when modeled by a uniform foreground screen. 

First, we consider a situation where the optical depth 
is not the same to all stellar particles, but the variations 
are uncorrelated with the intrinsic properties of the stel- 
lar particles. Such a scenario is by construction the case 
at the start of the simulation. For each stellar parti- 
cle individually the ratio of total to selective absorption, 
~ E{B-v) ~ 4-05, was taken from Calzctti et al. 
(2000). Since less extincted regions are also less red- 
dened and have a larger weight in the integrated SED, 
the effective extinction AveS = VAtt — Vint and effective 
reddening EiB - V%s = (b - V)Att - {B - F)i„t of the 
galaxy as a whole will not be related by the same factor 
4.05 as for the individual particles. Instead, the overall 
reddening for a given Ay will be smaller than predicted 
by Calzetti (i.e., the extinction is grayer). This is illus- 
trated in Fig. 10 where the dashed line represents the 
Av = 4.05 X E{B - V) scaling by Calzetti et al. (2000) 
and the dot-dashed curve represents a series of toy mod- 
els for (lifferc;nt Av^nmx with a uniform distribution of Ay 
values between and ^y,max to stellar particles that all 
emit at identical intrinsic luminosities. The effective vi- 
sual extinction for such a model consisting of N stellar 
particles emitting with identical luminosities Lj is calcu- 
lated as follows: 
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With AB,eS calculated in the same way as Av,eS in Eq. 
11, and using the relation between Ab and Ay to indi- 
vidual stellar particles from Eq. 7, the toy model gives 
an effective reddening of: 



E{B -VU = 2.5log 



0.4Av 



10 



(12) 

It is clear that the relation between total and selective 
absorption for the scries of toy models always lies below 
the Calzetti relation, and more so for the toy models with 
the largest Ay^eS (i-c, largest Av,max)- The reason is 
that stellar particles with little extinction are relatively 
dominant in the integrated observed SED. The binned 
distribution in Fig. 10, presenting total versus selective 
absorption for all considered simulation snapshots and 
viewing angles, is roughly encompassed by the Calzetti 
relation (Eq. 7) on the one hand and our series of toy 
models (Eq. 11 and Eq. 12) on the other hand. 

Since the Calzetti et al. (2000) attenuation law was 
derived empirically for galaxies as a whole, it is arguably 
not the appropriate law to apply to the individual stel- 
lar particles, i.e., the smallest stellar populations that our 
simulation can resolve, typically 10^ — 10^ Mq. We inves- 
tigated the changes in photometry when adopting a MW 
and SMC-like reddening curve by Pei (1992), which were 
derived in a more bottom-up fashion from the physics 
of interstellar dust grains. Again, we scaled the optical 
depth with the metallicity along the line of sight. For the 
SMC reddening curve, the resulting colors become red- 
der by up to 0.05, 0.1, and 0.2 mag in rest-frame B — V, 
U — V, and V — J respectively. The MW-like attenua- 
tion law is also less gray than Calzetti, thus producing 
slightly redder colors, though less so than for the SMC 
law. The effective extinction curve, expressed as ^(^B^y-j 
as a function of wavelength, of snapshots and viewing 
angles with large optical depths {Ay^gs > 1) is presented 
for different input attenuation laws in Fig. 11. In each 
case, we find the typical effective attenuation law to lie 
above (i.e., show less reddening for a given extinction 
than) the Calzetti et al. (2000) law. Even in the K- 
band, such an offset is still notable, and of similar size 
as would be expected for our toy model with Ay,cS = 2 
(corresponding to Ay^max = 7, i.e., many particles are 
heavily extincted). 



Not only does non-uniform extinction change the red- 
dening (^^) at a given Ay, it also affects the depen- 



dence of the reddening on wavelength C^^^-- 



For extinc- 



tion that is uncorrelated to the properties of the emitting 
sources, this gives the effective dust vector in the U — V 
versus V — J color-color diagram a shallower slope. I.e., 
for a given reddening inV— J, the reddening in [/ — F is 
smaller than predicted by the Calzetti et al. (2000) law. 
For example, for the simple toy models described in Eq. 
11 and Eq. 12, the slope ^'^•""--^v.eff ^j^^ effective red- 

dening vector decreases asymptotically with increasing 
Ay^max to half the slope of the Calzetti reddening vec- 



tor. The consequence of a different 



than Calzetti is 



clarified in Fig. 12. The solid curve represents the evo- 
lutionary track of a stellar population following a CSF 
history. The track starts 50 Myr after the onset of star 
formation and ends 2 Gyr later. Suppose different parts 
of a galaxy all contain a 1 Gyr old CSF population whose 
intrinsic location in color-color space is marked by the 
open gray box. A distribution of dust as described by 
the above mentioned toy model will redden the galaxy 
along the dotted line. Interpreting the observed colors 
{filled circle) as a CSF population attenuated by a fore- 
ground screen according to the Calzetti et al. (2000) law, 
will lead to a best-fit age {star symbol) that is too young 
and reddening that is too large. 

Since in our simulations the ages of the stellar parti- 
cles (that are each treated as SSPs) present at the start 
of the simulation were drawn randomly from a uniform 
distribution, the system has a CSF history in the earli- 
est snapshots without a correlation between the optical 
depth and intrinsic luminosities and colors of the stellar 
particles. Therefore, it comes as no surprise that, when 
looking at the attenuated stellar photometry in Fig. 13 
(for now all stars still set to solar metallicity), the cen- 
tral 68% interval in A log age^ reaches to more negative 
values (to -0.5 dex) during the earliest phases than was 
the case for the unattenuated photometry (Fig. 9). The 
estimated reddening is slightly larger than the true value, 
but nevertheless the use of Eq. 7 still causes an under- 
estimated Ay, as can be understood from Fig. 10. The 
systematic underestimate in age and Ay combined cause 
the evaluation of the stellar mass during the first snap- 
shots, when template mismatch due to the SFH is still 
negligible, to be too small by ~ 0.12 dex. 

After a few 100 Myr after the beginning of the simu- 
lation however. Fig. 13 reveals an improved recovery of 
the mass- weighted stellar age compared to that obtained 
by SED modeling of the intrinsic light (Fig. 9). Clearly, 
the assumption of a non-uniform dust distribution that is 
uncorrelated with the intrinsic properties of the emitting 
sources breaks down. 

4.3.2. Preferential extinction toward young stars 

Fig. 14 demonstrates the occurence of preferential ex- 
tinction toward young star forming regions in one of our 
simulations. The three panels indicate the binned dis- 
tribution of the metallicity-scaled hydrogen column den- 
sity measured along various lines-of-sight versus the age 
of the stellar particle to which the column density was 
computed for the 3 snapshots marked in Fig. 1(a). The 
vertical arrow indicates the start of the simulation. All 
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Fig. 11. — Effective extinction curves of simulated galaxies with Aygff > 1 for different input attenuation laws: (a) the Calzetti et al. 
(2000) law, MW-hke reddening from I'ei (1992), and (c) SMC-like reddening from Pei (1992). The black curve indicates the median over 
all snapshots and viewing angles with Ay ctt > 1- The hght gray polygon indicates the central 68% interval. The Calzetti, MW, and SMC 
atteimation laws are plotted in gray. In all cases, the effective extinction of simulated galaxies with large Ay^g is grayer than the Calzetti 
et al. (2000) law that is used in standard SED modeling. The offset is smallest when each stellar particle is attenuated according to the 
SMC-like law. 
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Fig. 12. — Rest-frame U — V versus V — J color-color diagram 
illustrating the effect of a non-uniform distribution of Ay values 
that is uncorrelated with the intrinsic properties of the emitting 
sources. The black curve indicates a CSF population witii age be- 
tween 50 Myr and 2 Gyr. Suppose an intrinsic population of age 
1 Gyr {box symbol) is reddened by such a dust distribution to the 
location in color-color space of the filled circle. Under the assump- 
tion of a uniform foreground screen of dust, the observed colors 
will then be traced back along the Calzetti et al. (2000) redden- 
ing vector {dashed black line), resulting in an artificially young age 
{star symbol). 



stellar ages older than this value (cut off for illustrational 
purposes) were set by hand as explained in §2.1. As we 
already pointed out in Fig. 1(f), the typical column den- 
sities are higher during the merger (panel b) than before 
(a) or after (c). Moreover, Fig. 14 shows that the ra- 
tio of column densities toward ongoing star formation 
over column densities toward older populations reac;hes 
a maximum during the merger (b) . Using sticky particle 
simulations of dusty starburst mergers, Bekki & Shioya 
(2001) found a similar age-dependent extinction, con- 
firming that this is a generic feature of merging systems 
and not determined by the method used to model dis- 



sipative processes. Poggianti & Wu (2000) inferred age- 
dependent extinction during a starburst to explain the 
nature of so-called e(a) galaxies: galaxies with [Oil] in 
emission and strong Balmer absorption lines, frequently 
associated with merger morphologies. 

From a physical perspec;tive, it is expected that during 
the merging process hydrodynamical and gravitational 
forces channel gas and dust to the central regions where 
it triggers a starburst. Once started, supernovae going of 
on a few 10^ yr timescale further increase the dust con- 
tent of the regions where newly formed stars reside. The 
fact that the distribution of younger (and thus intrinsi- 
cally bluer) stellar populations does not trace that of the 
older populations of stars and that it is intimately corre- 
lated with the dust distribution leads to an overestimate 
in age. In Fig. 12, a distribution of Ay values that was 
uncorrelated to the intrinsic properties of the emitters 
caused an effective reddening along a shallower slope in 
the U — V versus V — J diagram than for the Calzetti 
curve. In case of a distribution of Ay values that is cor- 
related with the intrinsic colors of the emitters (larger 
Ay to intrinsically bluer stellar populations), the slope 
of the effective reddening vector in the [/ — F versus V — J 
diagram is steeper than for the Calzetti et al. (2000) law. 
Since an observer will mistakenly model the galaxy with 
an intrinsically redder template, the reddening by dust 
E{B — V) will be underestimated. Although a given to- 
tal absorption corresponds to a stronger reddening in the 
presence of age-dependent extinction compared to uncor- 
related non-uniform extinction, the Calzetti et al. (2000) 
relation between E{B — V) and Ay given by Eq. 7 can 
still be considered as an upper limit. Therefore, the total 
absorption will be underestimated. This is illustrated in 
Fig. 10 where we plotted the binned distribution of true 
E{B — V) versus true Ay for all of our simulation snap- 
shots, viewed under a range of viewing angles. Finally, 
the derived stellar mass owes its more robust character 
to the compensating effects of systematic offsets in age 
and absorption. 

The effect of the larger extinction toward young stars 
will in practice be superposed on the effect of mis- 
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Fig. 13. — Effect of extinction. The difference between estimated and true (a) mass-weightod age, (b) stellar mass, (c) effective 
reddening, (d) effective visual extinction, and (e) star formation rate as a function of time since the merger. The SED modeling 
was performed on the attenuated stellar photometry with all stars set to solar metallicity. Style as in Fig. 7. Ages are still 
underestimated for the first 0.8 Gyr of the evolution, but to a lesser degree than estimates based on the intrinsic light. Added 
to the underestimated Av, this leads to a characterization of the stellar mass that is too low by 0.1 - 0.15 dex. The recovery of 
SFRs shows minima around first passage (~ —0.4 Gyr) and final coalescence (-0.1 to Gyr). 
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Fig. 14. — Distribution of hydrogen column densities, linearly 
scaled with the metallicity of the gas along different lines-of-sight to 
the stellar particles versus the age of the respective stellar particle. 
The relation between column density and stellar age is plotted for 
3 snapshots: before (2), during (9), and after (16) the merger (sec 
Fig. 1). The solid and dotted curves indicate the median and 68% 
interval of the distribution respectively. Darker intensity means a 
larger number of stars is present with that age. All stellar ages 
rightward of the arrow correspond to initial stars and were set by 
hand. The optical depth, which is proportional to the metallicity- 
scalcd gas density, is larger toward newly formed stars during the 
merger-triggered starburst. The signature of this age-dependent 
extinction weakens during more quiescent episodes of star forma- 
tion. 



match between template and true SFH, that prevents 
us from fully accounting for the difference between light- 
and mass-weighted stellar age (see §4.2). Fig. 15 illus- 
trates how an increased extinction toward young stars 
reduces the difference between the light-weighted and 
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Fig. 15. — Rest-frame attenuated and intrinsic V-band light- 
weighted age versus mass-weighted age for an initially 40% gas 
fraction simulation. The boxes mark the mean age weighted with 
the attenuated V-hand light. Darker intensities indicate a larger 
number of viewing angles. The solid and dotted curves mark the 
median and central 68% interval respectively. The dashed curve in- 
dicates the mean age weighted with the intrinsic V-band light (no 
attenuation) . The attenuated light- weighted age is a better approx- 
imation of the mass-weighted age than the intrinsic light- weighted 
age, increasingly so for younger stellar populations. Larger optical 
depths to young than to old stars are responsible for this effect. 



mass-weighted measure of age. We conclude that the 
SED modeling on galaxies with solar metallicity stars 

and dust distributed in between still underestimates the 
age during the merger-enhanced star-forming phases, but 
adding dust has improved our best guess to an overall 
median offset of -0.04 dex (compare Fig. 13(a) to Fig. 
9(a)). Similar conclusions were drawn by Bell & de Jong 
(2001) who examine the reddening and dimming effects 
of dust and its impact on estimating stellar mass-to-light 
ratios. At later times, the merger remnant is largely de- 
void of gas and dust, making accurate determinations of 
age and mass possible, as described already in §4.2. 

Finally, the evolution of AlogS'Fi? over time (Fig. 
13(c)) differs significantly from that obtained from mod- 
eling the intrinsic light (Fig. 9(c)). With the SED mod- 
eling now performed on attenuated photometry, we find 
underestimates in the SFR of -0.4 dex and -0.7 dex dur- 
ing first passage and final coalescence respectively. The 
agreement between recovered and true SFR improves 
during the first few 100 Myr after the starburst. Later, 
during the spheroid phase, large values of A log SFR 
merely reflect our inability to distinguish between low 
(a few 0.1 Mg/yr) and very low (~ 0.1 Mg/yr) SFRs. 
We therefore focus on the underestimates during the 
star-forming phases. In the presence of dust, there is a 
component of the merger-triggered starburst that is com- 
pletely enshrouded. This allows young stars to be hidden 
without revealing their presence by an increased redden- 
ing of the integrated SED. Since these heavily obscured 
star-forming regions are not accounted for by the best- 
fit model, the SFR during these phases is significantly 
underestimated . 



19 



2.0 


1 1 1 1 1 


1 1 i 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 

-i,-M,-,, -copulation 


1.5 


• Observed o/ 

it Recovered // 


1.0 






> 

1 




/ l-<y'' t300 




Z = 0.020 


0.5 




T/ Z = 0.008 - 






o 0.5 Gyr 


0.0 


o 1 Gyr '_ 
o 2 Gyr 


-0.5 


.... 1 


... 1 .... 1 .... 1 ... . 



-0.5 0.0 0.5 1.0 1.5 2.0 
V-J 

Fig. 16. — Rest-frame TJ — V versus V — J color-color diagram 
illustrating the eflfect of fitting solar metallicity templates to stellar 
populations of sub-solar metallicity. The black and gray curves rep- 
resent evolutionary tracks for an exponentially declining star for- 
mation history with e-folding time of 300 Myr for solar (Z = 0.02) 
and sub-solar (Z = 0.008) metallicity respectively, each starting at 
50 Myr. A stellar population with intrinsic colors indicated by the 
gray boxes will be reddened to the location in color-space marked 
by the filled circle in the case of (a) non-uniform age-independent 
extinction, (b) extinction by a uniform foreground screen, and (c) 
age-dependent extinction. In all three cases, the assumption of so- 
lar metallicity and Calzetti attenuation will lead to the conclusion 
that the stellar population formed its first stars 0.5 Gyr ago. This 
is an underestimate (a) or overestimate (b, c) respectively. The 
reddening is always underestimated. 



4.4. Impact of stellar metallicity 

In §4.2 and §4.3, we tested our SED modeling on syn- 
thetic photometry that was computed assuming a solar 

metallicity for all emitting sources. In reality, stars with 
a range of metallicities will be present, reflecting the level 
of enrichment in the gas at the epoch of their forma- 
tion. Before we repeat our analysis now setting the stel- 
lar metallicities to their appropriate value calculated by 
the GADGET-2 code, we anticipate the effect using the 
diagnostic U — V versus V ^ J color-color diagram in Fig. 
16. 

The tracks represent exponentially declining SFHs for 
metallicities of Z=0.008 (gray) and 0.02 (solar, in black). 
Both evolutionary tracks arc drawn from 50 Myr to 2 
Gyr after the onset of star formation. The classic age- 
metallicity degeneracy states that the optical broad-band 
colors of a young stellar population are nearly indistin- 
guishable from that of an older, more metal-poor pop- 
ulation (O'Connell 1986). For the raoo star formation 
history drawn here, this effect becomes only notable at 
later times: 2 Gyr after the onset of star formation the 
sub-solar metallicity track has the same U — V color as 
a solar metallicity population that started forming stars 
1.8 Gyr ago. On the one hand, the addition of dust will 
complicate the age-metallicity degeneracy. On the other 
hand, the addition of rest-frame NIR photometry helps 
to separate the evolutionary tracks for different metallic- 
ities (see Fig. 16). A galaxy whose attenuated light has 
colors marked by the filled circle may correspond with 
one of the intrinsic colors indicated by the gray boxes 
depending on the kind of extinction: (a) for non-uniform 



age-independent extinction, (b) for a foreground screen 
of dust, and (c) for age-dependent extinction. In case (a), 
the assumption of Calzetti attenuation and solar metal- 
licity in our SED modeling leads to a recovered evolu- 
tionary stage that is too young, marked with the star 
symbol on the solar metallicity track. In case (b) and 
(c), the same recovered evolutionary stage is too old. In 
all cases, the determination of the reddening will be too 
low, as will consequently be the case for the Ay and the 
stellar mass, and increasingly so for lower metallicities. 
Obviously, the effects described will again be superposed 
on the previously discussed effects of star formation his- 
tory and dust. It is also noteworthy that taking into 
account the enrichment by heavy elements reduces the 
effect of age-dependent extinction. Young stellar popu- 
lations are still intrinsically bluer than old populations, 
but to a lesser degree since they have formed at later 
times from gas that was more enriched. 

In our recovery analysis of stellar population proper- 
ties, we find that at metallicities of a quarter solar and 
below, the age is overestimated by to 0.5 dex (cen- 
tral 68% interval of Alogage,„). However, the under- 
estimate in reddening and therefore extinction for these 
low-metallicity galaxies is such that the mass estimate 
(which is dependent on both age and Ay) stays within 
±0.1 dex of its true value for 68% of the cases. We find 
no dependence of A log SFR on metallicity. 

4.5. Impact of AGN contribution 

Since the merger simulations described in this paper 
take into account the role of supermassive black holes on 

its environment (see e.g. Di Matteo et al. 2005; Springel 
et al. 2005b) , we now include its contribution to the inte- 
grated galaxy SED. We note that only during a timespan 
of the order of a Salpeter time, a few 10^ to 10^ year, the 
AGN emission amounts to a significant fraction of the 
stellar emission. Admittedly, the peak of AGN activity 
can be missed by the time sampling of our snapshots. 
Nevertheless, the current dataset provides a useful in- 
sight on its impact on the SED modeling. 

We illustrate the typical behavior in Fig. 17 showing 
the photometry computed at the time of merging when 
the accretion onto the supermassive black hole is max- 
imal. Here, the solid black curve represents the light 
received by an observer. We break down the attenuated 
SED in a stellar (dotted curve) and AGN contribution 
{dashed curve). Finally, the best- fit model (in this case 
an exponentially declining star formation history that 
started 0.8 Gyr ago) is plotted in gray. Although result- 
ing in a low Xreduced 

1, the SED modeling is mislead 
by a degeneracy between the stellar -|- AGN light and the 

stellar light of a younger population obscured by large 
columns of dust. The addition of AGN light, when ex- 
ceeding 10% of the total emission, adds another -0.1 to - 
0.15 dex to A log age^u, -1-0.05 to -hO.l mag to AE{B-V), 
and -1-0.3 to -1-0.5 mag to AAy. Since all observed light 
is treated as having a stellar origin, the star formation 
rates during the phase of and under viewing angles with 
significant AGN contribution to the observed SED are 
overestimated by an order of magnitude. These cases 
typically show a larger x^^duced (70% have Xreduced. > 5)- 
As expected, simulations with an initial gas fraction of 
80% have relatively more stars formed during final coa- 
lescence and more gas fed to the supermassive black hole 
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Fig. 17. — Attenuated spectrum of a simulated merger placed 
at z = 2.1 during the peak of AGN activity. The total attenuated 
light (black solid curve) is decomposed into a contribution from 
stars (dotted curve) and AGN (dashed line). An observer who sam- 
ples the total attenuated light with an identical set of broad-band 
filters as available for GOODS-CDFS and models the SED using 
stellar population synthesis only, will find as best-fitting model the 
spectrum in gray. Its age is too young by Alogagcu, = —0.4. The 
reddening E(B — V) and absorption Ay arc overestimated by 0.1 
and 0.4 mag respectively. The opposite sign of offsets in age and 
Av leads to a mass recovery that is only 0.05 dex below its true 
value. 



than lower gas fraction runs. Hence the bias due to a 
contribution of AGN light is more prominent for simula- 
tions with a higher gas fraction. 

4.6. Lessons for SED modeling 

How can the modeling of real high-redshift galaxies 
benefit from our analysis of merger simulations? In order 
to answer this question, we need to translate the mecha- 
nisms described above in terms of the physical properties 
of stars and dust into their combined effect as a function 
of observables. 

An observable that correlates crudely with time since 
the merger, the parameter function of which we 

quantified the quality of SED modeling, is specific star 
formation rate {sSFR = SFR/M). Generally, merger 
remnants are characterized by a low sSFR < 10~^-^ 
yr~^, whereas the sSFR is large {sSFR > 10^-^ yr~"^) 
during the merger-triggered starburst. However, very 
gas-rich isolated galaxies (the 'disk' phase of simulations 
with /gas,init = 80%) may reach similarly high sSFRs. 

Since the results of our SED modeling showed no trend 
with redshift for which the synthetic photometry was 
computed, we can also describe the quality of recover- 
ing stellar population properties in terms of rest-frame 
colors. Although not true observables, their computa- 
tion as described by Rudnick et al. (2003) suffers from 
only a minor template dependence compared to param- 
eters such as mass, age, dust content, and SFR. In Fig. 
18. we present the performance of our SED modeling as a 
function of location in the rest-frame U — V versus V — J 
color-color diagram. Simulations with initial gas frac- 
tions of 40% and 80% are shown separately. The SED 
modeling was based on the full (stellar -|- AGN) photom- 
etry in both cases. Downward triangles indicate a me- 



dian value of A = parameterrecovered — parametertrue 
that is negative for simulations with the respective col- 
ors. Upward triangles represent an overprediction of the 
true value. Lighter symbols are used for a better cor- 
respondence between the true and modeled parameter 
value. The different initial stellar ages and lower metal- 
licities for /gas = 80% runs explain why they extend to 
bluer U — V colors. At blue optical colors {U — V < 1), 
the attenuation is seriously underestimated. An overes- 
timate of the age by a similar order of magnitude (at 
U — V < 0.65 mostly due to the lower metallicity and 
at 0.65 < U — V < 1 due to age-dependent extinction) 
leads to a relatively robust estimate of the mass. 

At red optical colors {U — V > 1), sources with rela- 
tively blue V — J colors are generally better modeled than 
those with the reddest V — J colors in our sample. In the 
1 < U~V < 1.8 color regime, the effects of dust attenua- 
tion play an important role. Here we find objects that are 
heavily extincted during the merger-triggered starburst, 
but also disk galaxies seen edge-on during the earliest 
phases of the simulation. All sources with F — J > 1.5 in 
our sample belong the latter category. Since we did not 
impose an age gradient for the initial stars, the dust dis- 
tribution in thoses cases is non-uniform but uncorrelated 
with intrinsic color of the stellar particles. This explains 
why they have the strongest overestimate in E{B — V) 
(see §4.4). The upper part of the color-color distribution 
{U — V > 1.8) is where each galaxy ends up by aging 
without further infiow of gas. In such a system, the Ay 
values are modest to low and the lack of a template ex- 
actly matching the SFH is less problematic as the epoch 
of major star formation lies further back in time. 

5. RESULTS FROM SED MODELING WITH FREE 
REDSHIFT 

In practice, complete spectroscopic surveys of mass- 
limited samples of high-redshift galaxies are rare. Con- 
sequently, we often are not able to fix the redshift in the 
SED modeling procedure to its exact true value. Over 
the past few years, several codes have improved on es- 
timating the redshift based on the broad-band SED, by 
experimenting with different template sets and fitting al- 
gorithms. In this section, we will use the new photomet- 
ric redshift code EAZY (Brammer et al. in preparation) 
to establish the quality of photometric redshift (zphot) es- 
timates derived from our synthetic photometry and an- 
alyze the impact of 2;phot uncertainties on the derived 
stellar population properties. 

5.1. The photometric redshift code EAZY 

Here, we summarize the main characteristics of the 
EAZY photometric redshift code. A full description of 
the algorithm and template sets will be presented by 
Brammer et al. (in preparation). We test our ability 
to recover the redshifts using a template set based on a 
set of ~ 3000 PEGASE 2.0 (Fioc & Rocca-Volmerange 
1997) templates. The large number of templates was re- 
duced to 5 principal components necessary to span the 
colors of galaxies in the semi-analytic model by De Lu- 
cia & Blaizot (2007). In addition to these 5 templates, 
a dusty template with an age of 50 Myr and Ay =2.75 
accounts for the existence of dustier (and thus redder) 
galaxies than present in the semi-analytic model. 

The code fits a non- negative superposition of templates 
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Fig. 18. — Median quality of recovered stellar population properties for simulations with a gas fraction of 40% and 80% as a function of 
rest-frame optical and optical-to-NIR color. Downward triangles indicate underestimates with respect to the true value. Upward triangles 
mark overpredictions. Galaxies in regions with white triangles are characterized most accurately. Simulations with /gas = 80% reach 
bluer colors in {U — y)rest than those with /gas = 40% since their initial stars were set to younger ages (sec §2.1). At red optical colors 
[{U — y)rest > 1], galaxies with relatively blue (V — J)rest colors arc better recovered than those at the red (V — J)rcst end. The former 
are older, less obscured systems, while the latter have a young and dusty nature. At blue optical colors [(U — V')rest < 1], large systematics 
in the determination of age and Ay occur. Their opposite signs cancel out in the derivation of stellar mass. 



to the i3-to-8 /im photometry, using a template error 
function that effectively downweights the rest-frame UV 
and NIR portion of the templates in the fit. The value 
of the rcdshift marginalized over the total redshift prob- 
ability distribution is adopted as best Zphot estimate. 

EAZY allows for the use of a magnitude prior function, 
constructed from observed or simulated number counts 
as a function of apparent magnitude and redshift. How- 
ever, since we shift each simulation over the entire rcd- 
shift range 1.5 < z < 2.9, our methodology does not 
allow to test this feature. 

5.2. Recovering redshifts and stellar population 
properties from broad-band photometry 

5.2.1. Recovering redshifts 

In Fig. 19, we compare the photometric redshifts ob- 
tained with EAZY with the input redshifts for which the 
synthetic broad-band SEDs were computed. The mea- 
sure commonly used to quantify the photometric redshift 
quality is Az/{l + z). Its normalized median absolute de- 
viation is 0.031. We find no dependence of Az/{1 + z) 
on redshift, time since (or before) the merger, or mass- 
weighted age. A slight trend with Ay was found in the 
sense that the median Az/(l + z), which varies between 
-0.01 and 0.01 for < Ay < 2.5, increases to a few 
0.01 for simulated galaxies with Ay > 2.5. However, 
this systematic offset is still smaller than the scatter in 
Az/(1 -|- z) for these highly obscured galaxies. The per- 
formance of EAZY is very good and competitive with 
that of other codes presented in the literature. We note 
that the high ^phot quality owes greatly to the dense sam- 
pling of the rest-frame UV-to-NIR SED by our synthetic 
photometry (filter set and depths were adopted from the 
GOODS-Soiith survey). Leaving out 1 of the 11 filters, 
especially when located near a spectral break at the con- 
sidered redshift, may boost the uncertainty in ^phot- E.g., 
running EAZY on the synthetic SEDs without iJ-band 
results in a scatter over the 1.5 < z < 2.9 range that is 
35% larger: £7nmad(A2;/(1 + z)) = 0.042. At 2; ~ 2.7, 
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Fig. 19. — Comparison of photometric redshifts by EAZY versus 
true redshift. The correspondence is good (crNMAD(A2:/(l -I- 2:)) = 
0.031), without significant systematic offsets. 



where the Balmer/4000A break has shifted in between J 
and iJ, we find a median offset and scatter in A2;/(l -|- z) 
of 0.050 (zphot being larger than ^spec) and 0.071 respec- 
tively when omitting the H-band. 

This exercise offers a valuable complementary test to 
the empirical comparison with spectroscopic samples of 
high-redshift galaxies. The latter are direct measure- 
ments and therefore insensitive to our knowledge of stel- 
lar tracks and population synthesis. On the other hand, 
spectroscopic samples of high-redshift galaxies often suf- 
fer from selection biases, especially against galaxies lack- 
ing emission lines. 

5.2.2. Impact o/zphot uncertainties 
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Having quantified the quality of photometric redshifts, 
wc now repeat the SED modeling fixing the redshift to 
its best-fit value. The same mechanisms as discussed in 
§4 are still influencing the recovery of physical param- 
eters. Given the partially random nature of the 2;phot 
uncertainties, it comes as no surprise that the central 
68% interval broadens with respect to the SED modeling 
at fixed redshift. Averaged over time, the broadening 
is limited to 1%, 7%, 5%, 9%, and 1% in Alogage^, 
A log M, AE{B - V)es, AAv,eS, and A log SFR respec- 
tively, owing to the small scatter in 2;true versus Zphot- 

Any systematic part of the Zphot uncertainty translates 
into additional systematic offsets in the stellar popula- 
tion properties. Qualitatively, when a source is mistak- 
enly placed at higher redshift, a larger mass estimate 
and lower dust reddening are required to match the ob- 
served SED. Although we established that systematic 
offsets in Az/(1 + z) are small (§5.2.1), such an effect 
is noticeable for the most obscured {Ay > 2.5) galax- 
ies. The slight overprediction of their redshifts trans- 
lates into an additional offset in the estimated mass of 
logMfixz — logMphotz ~ —0.05 dex, in the estimated red- 
dening of t/)^^ -£^(S-F)photz 0.018 mag, and 
in the estimated extinction of Ay^x-^ — ^y.photz ~ 0.11 
mag. We note that, for an excellent photometric redshift 
accuracy as presented here, these effects are up to an 
order of magnitude smaller than the systematic offsets 
from the true stellar population properties as introduced 
by our poor knowledge of the star formation history, dust 
content and distribution, and metallicity. Furthermore, 
the impact of Zphot uncertainties discussed here is spe- 
cific for the particular template set used. Other template 
sets may introduce different redshift biases and therefore 
propagate into different biases when characterizing stel- 
lar populations. 

6. SUMMARY 

We analyzed the performance of a simple SED model- 
ing procedure applied to synthetic optical-to-NIR broad- 
band SEDs of merger simulations placed at redshifts 
z=1.5, 1.7, 2.9. The synthetic SEDs span a large 
range of colors and shapes, and arc significantly affected 
by dust attenuation. For all simulated galaxies, we find 
counterparts in the real high-redshift universe with simi- 
lar SED types. First, we modeled the SEDs assuming the 
redshift was known. The masses, ages, E{B — V), and 
Ay of simulated ellipticals are very well rcMsroduced, with 
an average value of AlogM = — 0.021q °j, Aloga^e^, = 
-Omto ll ^E{B -V) = -0.03j:°:J^, and AAy = 
— 0.291q;3q. Here the errors indicate the central 68% in- 
terval of the distribution of Aparameter values of all 
the simulations (different masses, gas fractions, viewing 
angles) in the spheroid regime. The simulated ellipticals 
are correctly classified as having low SFRs (A log SFR — 
— 0.23^0.47)- earlier, actively star-forming, phases, the 
scatter in recovered stellar population properties with 
respect to the true value increases, and larger system- 
atic underestimates of age, mass, extinction, and SFR 
occur. This is particularly the case for the simulation 
snapshots of phases with merger-triggered star forma- 
tion, where we find the following offsets and scatter 
(averaged over the merger regime indicated in Figure 
1): AlogM = -0A3l:°o\l Alogage^ = -0.12t°i°, 



AE{B -V) = -Omto oh AAy = -OM-toit and 
A log S-Fi? = -OAAto lv The SED modeling per- 
forms better on regular star-forming disks than on galax- 
ies during the merging event. Compared to spheroids 
however, the results of the SED modeling on disks 
show a larger scatter and larger systematic underesti- 
mates: AlogM = ~Omto% Alogage^j = -|-0.03tg;4^, 
AEiB -V) = -Omtom, ^Ay = -0.35lg:^|, and 

Alog S FR = -0.22toil. 

By adding the effects of dust attenuation, metallicity 
variations and AGN step by step to the basic intrinsic 
photometry, we were able to disentangle the different 
mechanisms at play and their impact on the estimation 
of the galaxy's stellar population properties. The quali- 
tative impact on the SED modeling results by different 
aspects of the galaxy content is summarized in Table 1. 

A mismatch between the real SFH and the allowed 
template SFHs leads to an inability to account for the dif- 
ference between light- weighted and mass-weighted prop- 
erties such as stellar age and mass. If the optical depth 
toward intrinsically bluer emitting sources is larger than 
to intrinsically redder stellar populations, the nett effect 
of the age and mass imdercstimatc due to mismatch be- 
tween true and template SFH is less severe. We find 
proof of such an increased extinction toward younger 
stars during the merger-triggered starburst. Applying 
the Calzetti et al. (2000) reddening law toward each stel- 
lar particle, wc find that the overall reddening for a given 
Ay is less than predicted by the Calzetti et al. (2000) 
law, particularly when the optical depth is uncorrelated 
with the intrinsic colors of the sources it is hiding. In the 
latter case, the effective dust vector has a shallower slope 
in the U — V versus V — J color-color diagram than the 
Calzetti et al. (2000) vector. In the case of larger optical 
depths toward young (blue) stellar populations, there is 
relatively more reddening in [/ — V for a given reddening 
inV — J. Applying a MW or SMC-like attenuation law 
to the individual stellar particles in the simulation in- 
creases the reddening, but the effective extinction is still 
grayer than the Calzetti et al. (2000) law. 

All other properties remaining the same, the effect of 
applying our SED modeling to stellar populations with 
sub-solar mctallicitics is that one would undcrpredict the 
reddening. For the young ages where such sub-solar 
metallicities could be expected, interpreting the light as 
coming from a solar-metallicity population will lead to 
an overestimate of the stellar age. 

Finally, our SED modeling is based on purely stellar 
emission. During the brief period when the AGN contri- 
bution is significant, the addition of its light will make 
the galaxy look younger and dustier, with a larger SFR. 

Although some of the biases described above could 
be anticipated from simpler models (see, e.g.. Fig. 12 
and Fig. 16), the simulations used in this paper prove 
very valuable in illustrating and quantifying the effects 
at play. Moreover, by following the formation of stars 
and distribution of gas and dust on a physically moti- 
vated basis, they provide insights that could not be ob- 
tained from simple toy models, such as the preferential 
extinction of young stars during the merger and its im- 
pact on recovering stellar population properties by SED 
modeling. 

We next repeated our analysis adopting the best-fit 
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Table 1. Qualitative summary of systematic trends in recovering stel- 
lar POPULATION PROPERTIES. 



Effect on age on E{B — V) on Ay 


on mass on SFR 


SFH mismatcli (higher SFR than in the past) - + 


+ 


Non-uniform extinction (uncorrelated with intrinsic colors of emitters) - + 




Age-dependent extinction (more extinction toward young stars) + - 


_a 


Z<Zq + - - 


+^ 


AGN - -t- + 





^The overestimated age and underestimated Ay compete in the determination of stellar mass. Since in pra<;tice we find this effect to be 
most outspoken in phases where the SFH mismatch is largest, the stellar mass will eflfectively be underestimated. 

''Here too, the systematic offsets in age and Ay have an opposite sign. From comparison of mass estimates based on the full attenuated 
photometry and its equivalent with all stars fixed to Zq , we find the effective A log M is positive. 

•^Similaj: to (b), we compared the results from SED modeling with/out AGN contribution and find that in the median the addition of 
AGN light increases the mass estimate slightly. 



photoinetric redshift estimate. Using the photoniet- 
ric redshift code EAZY (Brammer et al. in prepara- 
tion), we obtain a median normahzed absolute deviation 
cnmad(Az/(1 + z)) ~ 0.031. The random uncertainty in 
-2phot boosts the scatter in the quality measures A log M, 
Alog agen,, AE{B - V), AAy, and A\og SFR by less 
than 10%. A slight dependence of Az/{1 + z) on ex- 
tinction propagates into an additional systematic differ- 
ence in the estimated stellar mass, on the 10% level for 



galaxies with Ay > 2.5. Offsets in reddening and visual 
extinction Ay are anti-correlatcd with Az/{1 + z). 
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